
* Store Bootstrapped cluster SE results

* ---------------------------------------------	*
* 			Preparatory data					*
* ---------------------------------------------	*

* Create an empty storage .dta file to record results
clear
set obs 1
gen run=. 
save ${monte_data_dir}/storage_runs_cluster.dta, replace


* Prepare regression data to merge with bootstrap data
use BvDIDnumber year industry ///
lswMPm_1995_a hswMPm_1995_a vaempMPm_1995_a lngdpgap_1995_a gdppcMPm_1995_a ///
lswMPm_shr4_foreign_1995_a hswMPm_shr4_foreign_1995_a vaempMPm_shr4_foreign_1995_a lngdpgap_shr_foreign_1995_a gdppcMPm_shr4_foreign_1995_a ///
k${depvar}_${ttt} k${depvar}_${ttt}0 kNOT_${depvar}_${ttt} kNOT_${depvar}_${ttt}0 ///
spill${depvar}${ttt}_1995_a spillN${depvar}${ttt}_1995_a spill${depvar}${ttt}_1995_a0 spillN${depvar}${ttt}_1995_a0 ///
${depvar}_${ttt} ///
country_shr_1995 missing_weights_1995 maxweight_1995 missing_spill_weights_1995 ///
using ${final_dir}/regression_dataset${weight_window}${tf}.dta, clear

clonevar LSW = lswMPm_1995_a
clonevar HSW = hswMPm_1995_a
clonevar VAEMP = vaempMPm_1995_a
clonevar GDPGAP = lngdpgap_1995_a
clonevar GDPPC = gdppcMPm_1995_a

clonevar GDPPC_foreign = gdppcMPm_shr4_foreign_1995_a
clonevar VAEMP_foreign = vaempMPm_shr4_foreign_1995_a
clonevar GDPgap_foreign = lngdpgap_shr_foreign_1995_a
clonevar LSW_foreign = lswMPm_shr4_foreign_1995_a
clonevar HSW_foreign = hswMPm_shr4_foreign_1995_a

clonevar spilloversown = spill${depvar}${ttt}_1995_a
clonevar spilloversother = spillN${depvar}${ttt}_1995_a
clonevar spilloversownzero = spill${depvar}${ttt}_1995_a0
clonevar spilloversotherzero = spillN${depvar}${ttt}_1995_a0
clonevar stockown = k${depvar}_${ttt} 
clonevar stockownzero = k${depvar}_${ttt}0  
clonevar stockother = kNOT_${depvar}_${ttt} 
clonevar stockotherzero = kNOT_${depvar}_${ttt}0
drop lswMPm_1995_a hswMPm_1995_a vaempMPm_1995_a lngdpgap_1995_a gdppcMPm_1995_a ///
lswMPm_shr4_foreign_1995_a hswMPm_shr4_foreign_1995_a vaempMPm_shr4_foreign_1995_a lngdpgap_shr_foreign_1995_a gdppcMPm_shr4_foreign_1995_a ///
spill${depvar}${ttt}_1995_a spillN${depvar}${ttt}_1995_a spill${depvar}${ttt}_1995_a0 spillN${depvar}${ttt}_1995_a0 ///
k${depvar}_${ttt} k${depvar}_${ttt}0 kNOT_${depvar}_${ttt} kNOT_${depvar}_${ttt}0
mmerge BvD using ${final_dir}/bvd_list_regfirms_auto95.dta, unmatched(none) 
drop _m
keep if year >= 1995 & year <= 2011
sort BvD year
save ${monte_data_dir}/bvd_year_firmdata_clusterbootstrap.dta, replace

* ---------------------------------------------	*
* 			Run Boostrapping					*
* ---------------------------------------------	*

forvalue run = 1 / `monte_runs' {

	* Bootstrap data
	* ----------------

	di "Cluster: which run are we doing?"
	di "=== Cluster bootstrapping: `run'/`monte_runs' ==="

	di "cluster run `run' - bootstrap data"
	qui { 
		* Bootstrap the home countries
		use BvDIDnumber country_shr_1995 using ${monte_data_dir}/bvd_year_firmdata_clusterbootstrap.dta
		gsort BvDIDnumber -country_shr_1995	
		bysort BvDIDnumber country_shr_1995: keep if _n==1
		
		distinct country_shr_1995
		di "`r(ndistinct)' home countries"
		* the actual sampling (with replacement of course)
		bsample `r(ndistinct)', cluster(country_shr_1995) idcluster(country_id_new)
		gen lse_id_new = _n
		
		expand 17
		bys lse_id_new: gen year = 1994 + _n
		drop country_shr_1995


		* m:1 merge the bootstrap data with the regression data
		mmerge BvD year using ${monte_data_dir}/bvd_year_firmdata_clusterbootstrap.dta, unmatched(none)
		drop _m
		drop country_shr_1995
		ren lse_id_new lse_id
		ren country_id_new country_shr_1995

		* Prepare the regression indicators
		bys lse_id : egen _total_${depvar}_${ttt}_1995 = sum(${depvar}_${ttt}) if year>=1995+2 & year <= 2009+2
		bys lse_id : egen total_${depvar}_${ttt}_1995 = max(_total_${depvar}_${ttt}_1995)
		tab year, gen(YY_)
		egen yearctry = group(year country_shr_1995) if year <= 2009
		egen yearindustry = group(year industry) if year <= 2009
		xtset lse_id year
		sort lse_id year
	}

	* Run regressions:
	* ----------------

	* 1) LSW HSW GDPGAP | F+IY
	di "cluster run `run'- column 1"
	qui ppmlhdfe F2.${depvar}_${ttt} LSW HSW GDPGAP stockown stockownzero stockother stockotherzero spilloversown spilloversownzero spilloversother spilloversotherzero ///
		if missing_weights_1995==0 & maxweight_1995 < 1 & total_${depvar}_${ttt}_1995>0, absorb(lse_id yearindustry) vce(cluster country_shr_1995)
	local rc1=_rc
	local se_LSW_2 = _se[LSW]
	mat v = e(V)
	local stdLSW_1=v[1,1]^(1/2)

	di `se_LSW_2'
	di `stdLSW_1'

	if _rc==0 {
		mat b = e(b)
		mat v = e(V)
		local bLSW_1 = b[1,1]
		local stdLSW_1=v[1,1]^(1/2)
	}

	* 2) LSW HSW GDPGAP VAEMP | F+IY
	di "cluster run `run'- column 2"
	qui ppmlhdfe F2.${depvar}_${ttt} LSW HSW GDPGAP VAEMP stockown stockownzero stockother stockotherzero spilloversown spilloversownzero spilloversother spilloversotherzero ///
		if missing_weights_1995==0 & maxweight_1995 < 1 & total_${depvar}_${ttt}_1995>0, absorb(lse_id yearindustry) vce(cluster lse_id)
	local rc2=_rc
	if _rc==0 {
		mat b = e(b)
		mat v = e(V)
		local bLSW_2 = b[1,1]
		local stdLSW_2 = v[1,1]^(1/2)
	}

	* 3) GDPGAP + GDPPC | F+IY
	di "cluster run `run'- column 3"
	qui ppmlhdfe F2.${depvar}_${ttt} LSW HSW GDPGAP GDPPC stockown stockownzero stockother stockotherzero spilloversown spilloversownzero spilloversother spilloversotherzero ///
		if missing_weights_1995==0 & maxweight_1995 < 1 & total_${depvar}_${ttt}_1995>0, absorb(lse_id yearindustry) vce(cluster lse_id)
	local rc3=_rc
	if _rc==0 {
		mat b = e(b)
		mat v = e(V)
		local bLSW_3 = b[1,1]
		local stdLSW_3 = v[1,1]^(1/2)
	}

	* 4) GDPGAP | F+IY+CY
	di "cluster run `run'- column 4"
	qui ppmlhdfe F2.${depvar}_${ttt} LSW HSW GDPGAP stockown stockownzero stockother stockotherzero spilloversown spilloversownzero spilloversother spilloversotherzero ///
		if missing_weights_1995==0 & maxweight_1995 < 1 & total_${depvar}_${ttt}_1995>0, absorb(lse_id yearindustry yearctry) vce(cluster lse_id)
	local rc4=_rc
	if _rc==0 {
		mat b = e(b)
		mat v = e(V)
		local bLSW_4 = b[1,1]
		local stdLSW_4 = v[1,1]^(1/2)
	}	

	* 5) GDPGAP + VAEMP | F+IY+CY
	di "cluster run `run'- column 5"
	qui ppmlhdfe F2.${depvar}_${ttt} LSW HSW GDPGAP VAEMP stockown stockownzero stockother stockotherzero spilloversown spilloversownzero spilloversother spilloversotherzero ///
		if missing_weights_1995==0 & maxweight_1995 < 1 & total_${depvar}_${ttt}_1995>0, absorb(lse_id yearindustry yearctry) vce(cluster lse_id)
	local rc5=_rc
	if _rc==0 {
		mat b = e(b)
		mat v = e(V)
		local bLSW_5 = b[1,1]
		local stdLSW_5 = v[1,1]^(1/2)
	}	

	* 6) GDPGAP + GDPPC | F+IY+CY
	di "cluster run `run'- column 6"
	qui ppmlhdfe F2.${depvar}_${ttt} LSW HSW GDPGAP GDPPC stockown stockownzero stockother stockotherzero spilloversown spilloversownzero spilloversother spilloversotherzero ///
		if missing_weights_1995==0 & maxweight_1995 < 1 & total_${depvar}_${ttt}_1995>0, absorb(lse_id yearindustry yearctry) vce(cluster lse_id)
	local rc6=_rc
	if _rc==0 {
		mat b = e(b)
		mat v = e(V)
		local bLSW_6 = b[1,1]
		local stdLSW_6 = v[1,1]^(1/2)
	}

	* 7) GDPGAP | F+IY+CY | foreign
	di "cluster run `run'- column 7"
	qui ppmlhdfe F2.${depvar}_${ttt} LSW_foreign HSW_foreign GDPgap_foreign stockown stockownzero stockother stockotherzero spilloversown spilloversownzero spilloversother spilloversotherzero ///
		if missing_weights_1995==0 & missing_spill_weights_1995 == 0 & maxweight_1995 < 1 & total_${depvar}_${ttt}_1995>0,absorb(lse_id yearindustry yearctry) vce(cluster lse_id)
	local rc7=_rc
	if _rc==0 {
		mat b = e(b)
		mat v = e(V)
		local bLSW_7 = b[1,1]
		local stdLSW_7 = v[1,1]^(1/2)
	}
	
	* 8) GDPGAP + VAEMP | F+IY+CY | foreign
	di "cluster run `run'- column 8"
	qui ppmlhdfe F2.${depvar}_${ttt} LSW_foreign HSW_foreign GDPgap_foreign VAEMP_foreign stockown stockownzero stockother stockotherzero spilloversown spilloversownzero spilloversother spilloversotherzero ///
		if missing_weights_1995==0 & missing_spill_weights_1995 == 0 & maxweight_1995 < 1 & total_${depvar}_${ttt}_1995>0,absorb(lse_id yearindustry yearctry) vce(cluster lse_id)
	local rc8=_rc
	if _rc==0 {
		mat b = e(b)
		mat v = e(V)
		local bLSW_8 = b[1,1]
		local stdLSW_8 = v[1,1]^(1/2)
	}
	* 9) GDPGAP + GDPPC | F+IY+CY | foreign
	di "cluster run `run'- column 9"
	qui ppmlhdfe F2.${depvar}_${ttt} LSW_foreign HSW_foreign GDPgap_foreign GDPPC_foreign stockown stockownzero stockother stockotherzero spilloversown spilloversownzero spilloversother spilloversotherzero ///
		if missing_weights_1995==0 & missing_spill_weights_1995 == 0 & maxweight_1995 < 1 & total_${depvar}_${ttt}_1995>0,absorb(lse_id yearindustry yearctry) vce(cluster lse_id)
	local rc9=_rc
	if _rc==0 {
		mat b = e(b)
		mat v = e(V)
		local bLSW_9 = b[1,1]
		local stdLSW_9 = v[1,1]^(1/2)
	}

	* Incorporate original estimates 
	foreach var of numlist 1/9 {
		import delimited ${monte_data_dir}/original_estimates.csv,  clear varnames(1) case(preserve)
		sum estimateL if column==`var'
		local beta_orig_`var'=`r(mean)'
	}

	* Store run results:
	* -----------------

	clear 
	set obs 1
	gen run=`run'
	foreach var of numlist 1/9 {
		gen rc_`var'=`rc`var''
		if `rc`var''== 0 {
			gen beta_L_`var'=`bLSW_`var''
			gen tvar_`var'=(`bLSW_`var''-`beta_orig_`var'')/`stdLSW_`var''
			gen beta_L_orig_`var'=`beta_orig_`var''
		}
	}

	append using ${monte_data_dir}/storage_runs_cluster.dta
	drop if run==. 
	save ${monte_data_dir}/storage_runs_cluster.dta, replace	
}

